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THE FD-METHOD FOR SOLVING NONLINEAR KLEIN-GORDON 

EQUATION 

VOLODYMYR L. MAKAROV, DENYS V. DRAGUNOV, AND DMYTRO A. SEMBER 

Abstract. In the paper we present a functional-discrete method for solving the Goursat prob- 
04 ' lem for nonlinear Klein-Gordon equation. The sufficient conditions providing that the proposed 

^~^ ' method converges superexponentially are obtained. The results of numerical example presented 

^^ . in the paper are in good agreement with the theoretical conclusions. 

>^: 

! 1. Introduction 

in . 

^ ' It is well known that the Klein-Gordon equation (KGE) 

<[ (1.1) ^-^ + MW«.0)^*K.O 

r^ ' has extensive applications in modern physics and engineering. Particulary, it arises when study- 

^ • ing the scalar massive field in the de Sitter and anti-de Sitter spacetime [391 SO] , the propagation 

C ! of intense ultra-short optical pulses in low density dielectrics [13], the pionic atoms [27] et al. 

'—'I Furthermore, a partial case of the nonlinear KGE — the Sine-Gordon equation (SGE) — alone 

^ ■ has a great number of applications in physics. One encounters the SGE when studying the 

^ ■ propagation of a "slip" in an infinite chain of elastically bound atoms lying over a fixed lower 

chain of similar atoms [H], the magnetic fiux propagation in a large Josephson junction, the 
Q\ I domain wall dynamics in magnetic crystals |i2j et al. A nonlinear theory for strong interactions 

*^ ' has also been developed in which the SGE appears as a simplified classical model [301 EH]. In 

cn . geometry, the Goursat and Cauchy problems for the SGE are related to the existence of special 

nets on surfaces in £'^, which are called Chebyshev nets [32] . 

The whole range of the methods for solving the KGE can be conditionally divided into the 
two groups: the analytical methods and the discrete methods. The analytical methods allow us 
to express the exact solution to the equation through the elementary functions and convergent 
functional series. This methods a very useful for studying nonlinear physical phenomena, such 
as traveling waves and solitons [8J. Among the analytical methods there are the polynom,ial 
approxim,ation method ^ST[ i24j . the extended tanh method (see [25]), the sine-cosine method 
(see [3S]), the variational iteration method (see [31 S]), the homotopy methods (see [331 [D ESj 
and references therein), the infinite series methods and many other methods and techniques 
(see [361 1311 HH] and references therein). However, when the quahtative analysis of a solution 
is not the main target, then discrete methods can be useful as well. The methods of this 
group approximate the exact solution on a finite set of distinct points (i.e., on the mesh). The 
discrete methods for solving KGE take their origins mainly from the finite-difference methods 
(see [5l [291 M ^^^ references therein) and Runge-Kutta methods (see [281 C] and references 
therein) . 
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Apparently, the border line between the two groups of methods introduced above is very fuzzy 
and some methods can be considered as belonging to both groups at once. In the present paper 
we offer a numerical-analytical method which possesses the mein properties of both analytical 
and discrete methods simultaneously. This method (hereinaftr referenced to as the FD-method) 
is based on the FD-approach described in [121 EI] and takes its origins from the functional- 
discrete method for solving Sturm-Liouville problems (see [22| [23]). 

The paper is organized as follows. The Goursat problem for nonlinear KGE is introduced 
in Section |2l Section [3] is devoted to the description of the FD-method's algorithm for the 
given Goursat problem. In Section H] an important auxiliary statement about approximating 
properties of a hyperbolic differential equation with piece- wise constant argument is proved (see 
Theorem 14. ip . Theorem 14.11 plays a key role in the proof of Theorem 15.11 containing sufficient 
convergence conditions for the proposed FD- method (Section E]). A numerical example and 
conclusions are presented in Sections El and [7] respectively. 



2. Problem statement 
Let us consider KGE (II. ip in a slightly modified form 

(2.1) ^^^^ + Niu{x,y)) = f{x,y), 

oxoy 

which is more suitable for application of the FD-method. Equation (12. ip can be obtained from 
(II. ip via the transform of variables 

(2.2) t = x-y, i = x + y, 

assuming that 

u{x, y) =v{x-y,x + y), f{x, y) = ^{x-y,x + y). 

Since the FD-method is a numerical-analytical method (not mere analytical) it cannot be 
applied to equation (12. ip without an initial or boundary condition. In the paper we confine 
ourself to considering a Goursat problem (see [S]), supplementing equation (12. ip with the 
following boundary conditions: 

(2.3) u{x,0)=ij{x), u{0,y) = (l){y), ^(0) = 0(0). 
We assume that nonlinear function N{u) can be expressed in the form of 

oo 

n{u) = N{u)u, N{u) = Y^ UsU', z/, G R, Vti e R 

s=0 

and 

^(x) G C(i) (Di) n C (Di) , 0(y) G C(i) (D2) n C (D2) , fix, y) G CiD)E 

D = {{x,y):0<x<X,0<y<Y}, Di = (0;X), D2 = iO;Y) . 

Given assumptions imply that the solution u{x,y) G C^'^{D) fl C{D) to the Goursat problem 
(12. ip . (12. 3p exists and is unique (see [18j). 



Hereinafter a horizontal bar above a letter indicates the closure of the set denoted by the letter. 
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3. General description of the FD-method's algorithm for solving KGE 

According to the FD-method's algorithm for solving operator equations described in [12], the 
FD-method for solving Goursat problem fl2.ip . fl2.3p can be constructed in the following way. 

We approximate the exact solution u{x,y) to problem f l2.ip . f l2.3p by the function u[x,y) 
defined as the finite sum 

m 

(3.1) u{x,y) = ^ u{x,y), 

fc=0 

where m G N. In the rest part of the paper the function u{x, y) can be also referenced to as the 
FD-approximation of rank m. 

To define the functions u{x, y) we have to introduce a mesh 

X Y 

(3.2) Xi = hii, yj = /laj, hi = —-, /i2 = -— , 

iVi iV2 



zeO,iVi, jeO,N2, NuN2>l. 

For a while we assume that the positive integers A^i and N2 are chosen arbitrary, however, later 
it will be shown that decreasing the value of parameter h = ^/h^ + h^ (that is, increasing both 
A^i and A''2) we can increase the accuracy of the FD-method. 

As soon as mesh (13. 2p is fixed we can define function u[x, y) G C{D) as the solution to the 
nonlinear Goursat problem with piece-wise constant argument (hereinafter referenced to as the 
basic problem) 

(3.3) ^A^ + N(S{x,.„y,.,)) ^S{x,y) = f{x,y), \f{x,y) e P.„ 

(3.4) ^Six,0)=ij{x), u{0,y) = cP{y), ^'(O) = 0(0), V(x,2/),gD, 
where 



(3.5) Pij = {xi_i,Xi) X (?/j_i, yj) , i G 1, iVi, j G 1, A^2- 



H 



Once the basic problem (13. 3p . (13. 4p is solved, the functions u[x,y) G C{D), k E l,m can be 
found as the solutions to the following sequence of linear Goursat problems 

(fc) 

(3.6) ' ^ N{u{xi-i,yj-i)) u{x,y) + N'{u{x,y)) u{x,y) u{xi^i,yj^i) = 

Ak-s[N; u{xi-i,yj^i),... u [Xi-x,Vj-\)) u{x,y)- 

s=l 

~Ak{N; u{xi_i,yj_i),..., w (xi_i,?/j_i), O) u{x,y) + 



k—l 

I (0), {k-\-s) , 

Afe-i-^A^; M(xi_i,2/j_i),..., u {xi-i,yj-i))' 



s=0 

/ (0) (fc-l-sl -X 

-Aj^i_s[N;u{x,y),..., u {x,y)) 



u{x,y) =F{x,y), 



y{x,y)eP,,„ \/iel,N,, VjGl,A^2. 
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(k) (k) 

HO, y) = V(x, 0) = 0, Vx G [0, X] , Vy G [0, Y] . 



(3.7) 

Here v4„ (A^; vq, f i, . . . , Vn) denotes the Adoinian polynomial of n-tli order for the function A^(-) 
(see, for example, [37], [IE], [20]), which can be calculated by the formulas 



1 d"' ( 



n\ dr^' 



(3,i 



E 'v<'-.>(„„)^ 



ai— 02 



T=0 



ai+...-|-a,i=n 

ai>...>o„+i=0 

aiSNUJO} 



(«! -a2)! '" (ttn -a„+i)!' 



4. Approximating properties of the basic problem 

It is well known (see, for example, [6]) that problem (I3.3p . (13 ■4p possesses a unique solution, 
which can be represented in the following form: 

u[x,y) = R{x,yj-i,x,y) u{x,yj-i) + 

+R{xi^i,y,x,y) u{xi-i,y) - R{xi^i,yj-i,x,y) m(xj_i, 2/j-i)- 



(4.1) 



y 



d_ 



d 



R{^,yj^i,x,y) 



—R{xi^i,r],x,y) 



u{^,yj-i)d^- 



u{xi^i,r])dr]+ 



Vj-l 
X y 



R{i, V, X, y)f{^, v)dCdv, V(x, y) G Pi, 



3-> 



Xi-i Vj^i 



where 



(4.2) 



R{x, y- e, Ti) = Jo UANij{x-0{y-^y\ = oFi (1; -(x - Oiv - v)Ni,j) , 
—R{x, y- e, ?7) = oi^i (2; -(x - (){y - v)Nij) Ni^,{r, - y), 



d_ 
dy 



R{x,y;^,r]) = qFi (2; -{x - Oiv - v)Nr,,) iV,,,(e - x) 



N^,J 



,(0) 



and Jq, qFi denote the Bessel function of the first kind and the confluent hypergeometric function 
respectively (see [17]). Using integration by parts we can rewrite formula ( 14. ip as follows 



(4.3) 



u{x,y) =u{xi^i,y)+ I R{C,yj~i,x,y) 



' d (0) 



di- 



Xi-l 



d_ 
dr) 



R{xi_i,r],x,y) 



u{xi_i,r])dr]+ 



j/j-i 
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X y 

^i-l Uj-l 



Theorem 4.1. Suppose that u{x,y) and u[x,y) are the solutions to problems (12. ip . fl2.3p and 
(13. 3p . (13. 4p respectively. Then for the sufficiently small values of hi and /i2 there exists an 
independent on hi and /i2 constant k, such that 



(4.4) 



\u{x, y)— u{x, y)\\^ < hK, h = \jh\ + K. 



Proof of Theorem 4^1 ■ Let us consider the auxihary function 



z{x,y) = u{x,y)- u{x,y). 
It is easy to see that this function is continuous on D and satisfies the equation 



(4.5) 



d'^z(x,y) ^^,(0) ,, , , 

Q^Qy + ^( u{xi.i,y,-i))z{x, y) + 



+ 



M 



u{x,y) = 0, y{x,y)ePij, 



N{u{x, y)) - N{ u{xi_i,yj_i)) 

together with the boundary conditions 

z{x,0) = 0, z{0,y)=0, WxeDi,yeD2. 

To prove the theorem it is enough to find a positive real constant k, independent on h, such 
that 

(4.6) 



\z{x, 



P^. 



< hn. 



yiel.Ni, Vj G 1,A^2, where ||2;(x,|/) lip = max \z{x,y)\. 
For further convenience we have to introduce the notation 



(4.7) 



Uij=u{Xi_i,yj_i), Ui^s = u{Xi_i^y^_^), N^^. 



N'(l,) 



(4.8) 
(4.9) 



Li i = max 

{x,y)ePi.j 

9e[o,i] 



^'{uij -0{uij -u{x,y))) 



Zi,j = z{xi^i,yj^i), \\z\\-^j = \\z{x,y)\\p^^^ , 



(4.10) 



Ri,, = oFi(l; Nijhih^), B!ij = oFi{2- Ni,,hih2)Ni^j, 



(4-11) ||m|| = ||m(x,2/)||5,||^'|| = ||V''(2;)ll[o,x]'ll0ll = lkiy;ii[o,y]- 

It is worth to emphasize, that the principal role in the proof is assigned to the constants N^ 
and La defined in the following way 

(4.12) Na= max |A^(m)|, L„ = max |N'(m)|, 

u€{pi,p2) u€{pi,P2) 

Pi = min_ u{x, y) — a, p2 = max_ u{x, y) + a. 

{x,y)GD {x,y)£D 
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Here a denotes an arbitrary positive real number fixed throughout the proof. It is easy to see 
that according to the definition of Na f!4.12p we have the inequality \\N{u{{x,y))\\jj < N^- 
To prove Theorem 14. II we need the following auxiliary statement. 

Lemma 4.1. Suppose that u{x,y) and u[x,y) are the solutions to problems (12. ip . (12. 3p and 
(13. 3p . (13. 4p respectively. Then the following inequalities hold true 



kii,,<ikii,_i,(i + M2^y + 



i-i 



(4.13) 



[ s=l 

+/ii/i2-Rij {Lij + Nij) ||2;||jj_i + hih2hRij {Nij + Lij) A, 



+ 



(4.14) +h2R^Ah^J2[i'^Ns,J 



+ Ls.i) \\z 



s,JJ ll^llsj 



h {Nsj + L,j) A 



+ 



+hih2Rij {Lij + Nij) \\z\\._^ ■ + hih2hRij (Nij + Lij) A, 



for all {x, y) E Pij, i E 1, Ni,j E I, N2, wh 



ere 



(4.15) 



A 



Y(\\f\\ + N^\\u\ 



+ 



X(\\f\\ + N^\\u\ 



21 2 



loj 



z\\,, = 0, \/tEl,Ni,WjEl,N2. 



Proof of Lemma \4^ Unless otherwise stated we assume that (x, y) E Pij for some fixed positive 
integers i E l,Ni and j E 1, A^2- 

We begin with the proof of inequality (I4.13p . As it was mentioned above, function z{x, y) can 
be represented by virtue of the Reimann function in the following form: 

z{x,y) = z{x^^l,y) + 



(4.16) 



+ / R{C,yj-i,x,y) 






d^- 



y 



d 

—R{xi_i,r],x,y) 



z{xi^i,r])dr]+ 



%-i 



(0). 



N{u{^,r])) - N{u{x,.,,y,-i)) 



X y 

Equality (I4.16P yields the estimate (see notation (14. 7p - (14.101) ) 

dz{x, yj^i 



u{^,r])d^dr]. 



\z\\,. < \\z\\.^,. + hiRij 



dx 



hih2R'ij\\z\\ + 



(4.17) 






i,j 
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< 



/(O) X (0) -(0) , 



«j 



<iizii.„ . (i + /ii/i2/2y + /ii/2i„ 



+/ii/i2-Ri,j {Lij + Nij) \\u{x, y)- 



dx 



+ 



*j-i 



(0) 

u 



«.JI 



< 



*J 



< lkll,-ij(l + ^i^2i?y)+/iii?ij 



dz{x, Vj. 



dx 



+ 



i,j-'^ 



+hih2Rij (Lij + Nij) \\u{x, y) - Uij\\^^. + |zjj| 

Taking into account the obvious inequality l^jJ < |U|L,_i ,D we can summarize inequalities 
( I4.17P in the following way 

dz{x,yj^i] 



\4i,j<\\4^-iA^ + ^^^^Kj)+h^R^>J 



dx 



+ 



i,j-^ 



(4.18) 



+ /ii/i2-Rjj {Lij + Nij) \\z\\i,j^i + 
+hih2Rij {Lij + Nij) \\u{x, y) - UijWi^ 



Let us estimate the expressions 



dx 



i,j~^ 



and \\u{x,y) — u{xi_i,yj^i)\\. . arising in the 



right side of inequality f l4.18p . We start with the latter one. The mean value theorem provides 
us with the inequality 

\u{^,r]) -M(x,_i,yj_i)| < \u{x,y) ~u{xi^i,y)\ + 

+ \u{xi^i,y) - u{xi^i,yj^i)\ < 



(4.19) 



< \x — X. 



i-l\ 



du{x, y) 



dx 



+ \y-yj-i\ 



D 



du{x, y) 



dy 



D 



Furthermore, equation 02. ip yields us the equalities 



(4.20) 



du{x,y) 
dx 



ip'ix)+ / [f{x,r]) - N{u{x,ri))u{x,r])]dr], 



(4.21) 



du{x, y) 
dy 



<t>'{y)+ / [f{U)-N{u{U))n{U)]di- 



Combining (I4.19P with (I4.20p and (I4.2ip we get the estimate 

\u{x,y) -u{xi.i,yj^i)\ < 



(4.22) 



< h, (||V''(x)||[o,^] +Y[\\f{x,y)\\c + \\N{u{x,y))\\^ \\u{x,y)\\c] 
+h2[U'{y)%^y^+X[\\f{x,y)\\j, + \\N{u{x,y))\\a\\u{x,y)\\c]\< 



<A\hi + hi = Ah. 



The inequality \zij\ < \\z\\^_-^ ■ is valid as well and it will be used in the proof of inequality (|4.14p . 
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Now let US pass to the estimation of 



dz{x,yj^i) 
dx 



. Equation (14. 5p implies the equality 
dz{x,yj-i) _ 



dx 



_i ys 



(4.23) 



s=l. 



,(0) 



(0) 



~X^ / {N{ui^s)z{x,ri)+ N{u{x,ri)) - N{ui^s) u{x,ri)\dr] 



?/s-l 



Using notation (14.81) from (14.231) it is easy to get the estimate (Vx G [xi_i,Xj]) 
(4.24) 



dz{x, Vj-i 



dx 



<^Ni^s / \z{x,r])\dr]+ 

Vs-l 



i-1 



s=l. 



+ ^ / {Li^s + Ni^s){\u{x,r]) -Ui^s\ + \zi,s\} dr]. 



ys~i 



Combining inequalities (I4.22p and (I4.24p we obtain 
dz{x, Vj-i 



(4.25) 



dx 



s=l s=l 



Finally, using estimates (I4.22p and (I4.25P we get from (I4.18P the target inequality (I4.13p . 
The proof of inequality (I4.14p is mostly similar to the proof of inequality (I4.13p . However, to 
obtain inequality (I4.14p we have to use the formula 

z{x,y) = z{x,yj^i) + 



(4.26) 



+ 



d_ 



R{^,yj~i,x,y) 



z{i,yj-i)di- 



d_ 
drj 



z{xi^i,r]) 



(0) 



dr]+ 



y 
- / R{xi_i,r],x,y) 
yj-i 

X y 

+ Ri^^V,x,y) 

^i-i yj-i 
instead of formula (14. 16p . Formula (14.260 leads us to the inequality 

dz{xi-i,y) 



N{u{^,r]))-N{u {xi_,,yj_i)) 



u{^,v)dCdv- 



\z\ 



„-<lk||,,_l(l + /^l/^2i?U+^2i?,, 



dy 



+ 



«-iJ 



(4.27) 



+ /ii/i2-Ri,j {Lij + Nij) \\z\\^_^ ■ + 



+hih2Rij {Lij + Nij) \\u{x,y) - Uij\\^^. . 
Instead of inequality (I4.25P we have to use the following one (Wy G [yj-i,yj]) : 



(4.28) 



dz{xi^i,y) 
dy 



<hJ2 ('^^s,j + Ls.i) \\z\\^^^ + h^hA^{N,j + L 



S,J) 1 
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which can be obtained in a similar way. Finally, estimates fl4.22p . fl4.27p and (I4.28P lead us to 
inequality ( I4.14p . which was to be proved. 

The proof of Lemma 14.11 is completed. D 



To use the results of Lemma H?T] we have to make an assumption about the correlation between 
hi and /i2. It is precisely this assumption that determines which of the estimates, f l4.13p or fl4.14p . 
will be used in the further reasoning. Without loss of generality we assume thaio 



(4.29) 



hi<h2^ N2< 



X 



Taking into account estimate fl4.13p together with the definitions of constants iY^, and L^, 
(14.121) . we can conclude that inequalities 



,(0) 



(4.30) pi <Uk^i< p2, VA; G 0, iVi - 1, V/ G 0, j, ] < N^ 

imply the estimates 

\\z\\,,<{l + hh2K)\\z\\,_,,+ 



l-l 



(4.31) 

where 
(4.32) 



+ hiRa{2Na + L^)h2 J2 II^IL,s + hih2Ra{La + N^) |k||fc,;_i + 



s=l 



+hihRaA{Na + La){Y + h2), ke l,A^i, V/ G l,i + l. 



Ra = oFi (1; Nahl) , R'^ = oFi (2; N^hl) N^ 



However, generally speaking, conditions (I4.30p could not be satisfied for all/ G 1, A^2 unless 
some restriction on the value of ^2 is imposed. To find out this restriction we have to study 
some properties of the sequence of real numbers fiij Wi E 0, Ni, Wj G 0, N2 defined by formulas 

(4.33) fiij = afii^ij + 6/iij_i + c, /ioj = /ij,o = 

with 

a = 1 + /iiai(/i2) = 1 + /ii/i2-R^, 



(4.34) 



b = hMh2) = h RaY (2iV„ + La) + /i2^a {Na + La) 



c = hihci{h2) = hih RaA {Na + La) {Y + /l2) 



Lemma 4.2. Suppose that real numbers fiij Wi G 0,Ni, Vj G 0, A^2 c^^e defined by formulas 
(14331) with 

a = l + hiai, b = hibi, c = hihci 
and assumption (14.290 holds. Then 



(4.35) 



Hij < hXciexp({X + Y)bi + Xai\ Wi G l,A^i, j G 1,N2. 



^Under this assumption inequality (|4.14p is not useful for us and we need to use inequality (|4.13p instead. 
However, if we were assumed that hi > ft.2 we would be forced to use inequality (|4.14l) to prove the theorem. 
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Proof of Lemma \4.S\ Using the method of mathematical induction, it is not hard to prove the 
exphcit formula for calculation of fiij 



iii = 0,j el,N2 0i j = 0,iel,Ni, 

k=Op=0 

Using formula f l4.36p and assumption f l4.29p we get 

liij < fJ^NuN, = c}_^d' }_^ b = 

p=0 k=0 '^' 



p=o fc=o ■ V 1/ 

p=0 k=0 ^ ^^ 

p=0 fc=0 ^ / \ 1/ \ 1/ 

Afl-1 ^2-1 //-^ ^X , Nfc 



p=0 A:=0 

1 + — ai J exp ({X + Y)bi] < hXci exp Ux + y)foi + Xoi' 

The proof of Lemma 14.21 is completed. D 

To obtain the required restriction on the maximum value of h2 mentioned above, we have to 
consider the auxiliary function 

E{hi,h2) = ^/h[+hiciih2)Xexp(^{X + Y)bi{h2) + Xai{h2)y 
Taking into account assumption f l4.29p we arrive at the inequality 

E{huh2) <E{h2,h2) = £{h2). 
Function £{h2) is strictly increasing function on [0,+oo], lim £{h2) = +00, £^(0) = 0. This 
fact implies the existence and uniqueness of the constant Ha > 0, such that 

£{Ha) = a, and £{h2) < a,V/i2 G [0,Ha) ■ 
Henceforward we assume that 

(4.37) h2<Ha. 

It follows from Lemma [4.21 that inequality fl4.37p provides the estimates 

(4.38) fiij <hK <a, K = ci{Ha)X exp({X + Y)bi{Ha) + Xai{H^ 



iel,N^, je l,N2. 

Now we are in position to prove inequality (14. 6 p with k = K (see notation (I4.38P ) via the 
method of mathematical induction. We will use induction with respect to j. 
The base case, j=l. 
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Taking into account that u{x,0) =u{x,0) V x G [0,X] we arriveat the conclusion that 
inequahties fl4.30p are vahd for j = 0. This fact provides inequahties (I4.3ip for j = 1, which, 
using notations f l4.32p and f l4.34p . can be represented in the form of 



(4.39) \\z{x,y)l^^<\\zl_^^^a + c, yi,el,Ni. 

Taking into account formula f l4.33p it is easy to see that 



(4.40) \\z{x,y)l^,<fi,,i, V2Gl,iVi 

and, consequently, from (I4.38P it follows that 



(4.41) \\zix,y)l,<hK, VzGl,iVi 



Inequalities (14.41 p prove inequality (14. 6 p with k = K for j = 1 and for alH G 1, A^i. 
Induction step. 
Assume that inequality ( 14. 6p and auxiliary inequality 

(4.42) \\zix,y)l^.<fi,j 



are proved for all j G l,n, i G l,A^i, 1 < n < A^2- This assumption implies the validity 
of inequalities (l4.3Up for j = n and, consequently, we obtain inequalities (I4.3ip for j = n. 
Combining the auxiliary inequality (I4.42p together with obvious inequalities 



and inequalities (I4.3ip with j = n, we arrive at the following estimate for Zi^n+v 

n 

(4.43) + hiRa{2Na + L„)/l2 Yl ^*.« + hlh2Ra{La + Na)fli,n < 

s=l 

< a/ij_i^n+l + bfli^n + C = /ij,n+l- 



Estimate (I4.43P is valid for alH G 1, A^i and proves inequality (I4.42p for j = n + 1. Furthermore, 
taking into account estimates (14.380 and estimate (14.430 we immediately obtain the validity of 



inequality (14. 6 p with k = K, for j = n + 1 and for alH G 1, A'^i, which was to be demonstrated: 

\\^i^^y)\\i,n+i < ^^i,n+l <hK <a, V2 G 1, A^i. 

By the principle of mathematical induction it follows that inequality (14. 6 p is valid for all 
z G 1, A^i, j G 1, A'2. Thereby, the theorem is proved. D 

5. The FD-method for nonlinear Goursat problem: the convergence result 

In this section we intend to study the question of sufficient conditions providing convergence 
of the FD-method (O), (I32D, (O, (133!), (ESD, (ETD to the exact solution of the Goursat 
problem (12. ip . (12. 3p . In other words, given that the parameters hi and /i2 are sufficiently smaljj 
we will prove that 

m 

(5.1) lim y]|U(2:,2/)|L /5 < 00, 

k=0 



We assume that hi < h2 < I and inequality (49) holds true. 
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and 

(5.2) 

where 



(kl 



uix,y) = ^ u{x,y) 



l|/(^'?/)|li,S= 



max 



||/(2^>1/)| 



D' 



A;=0 



m ax 

i&l,Ni,jel,N2 



i^^(^'^)iik. + ii|^(^'^)iik 



for all f{x, y), such that /(x, y) E C{D) and /(x, y) G C^'"^(Pjj), i G 1, A^i,i G 1, A^2- To achieve 
that we will use the method of generating functions and the main part of this section is devoted 
to the derivation of an appropriate equation for a generating function. 



{k\ 



We begin with the estimation of m(x, y) 



D^ 



W, 



The piece-wise constant function u^{x^ y) , (x, y) E D defined as 



(^), 



W. 



u±{x,y) =u{xi_i,yj_i) , if {x,y) G [xi^i,Xi) x [y^.i/i/j), Vi G l,A^i,Vj G 1, A^2 
allows us to represent equations (13. 6p in the form which is valid for all (x,y) G D, i.e., 

d^ u(x,y) ,,/(o), ,x (fc), , 
+ A^'(m^(x,?/)) u{x,y)u±{x,y) = 



(5.3) 



fc-i 



(0) 



(fc-s)^ 



w. 



^Afc_4A^;w±(a;,?/),... w±(x,i/)) w(x,?/)+ 



s=l 

fc-1 

s=0 



^fc-i-s(A^,^±(a:,y),..., w± ix,y))- 



(0). 



-Afc-i--s(A^, w(a;,2/) 



,..., u {x,y)) 



u{x,y)- 



(0), , (fc-l)^ , ^^ (0)^ , (fc). 



-Afc(A^,M^(x,?/),..., M±(x,y),0) u{x,y) =F{x,y). 

As it was mentioned above (see previous section or [19j), the solution to the fc-th equation of 
system (13. 6p {k > 1) on Pj j can be represented in the following form: 



(5.4) 



u{x,y) =u{xi-i,y)+ I R{C,yj-i;x,y) 



d (k) 



d^- 



Xi-\ 



d_ 
di] 



R{xi.i,r]]x,y) 



Vj-i 



,.(0) (fc) 



X y 



-N'{u,,) 



u 



*j 



u{xi^i,r])dr]+ 



H 



R{^,V;x,y) u{^,r])d^dr]+ 



Xi-i Vj-i 
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(0). 



(0). 



As it follows from Theorem 14. ![ the function u(h,x,y) =u{x,y) tends uniformly to u{x,y) 
on Z) as /i —7- 0. Hence, taking into account the existence and uniqueness of the continuous on 
D solution u{x,y) to the Goursat problem (12.11) . (12. 3p . we can conclude that there exists an 
independent on hi and /i2 constant M„, such that 



(5.5) 



u{x,y) 



D 



<M„. 



The last fact provides the existence of the independent on hi and /12 constants M^, M'^, Mr, M'j^ > 
0, such that 



(5.6) 



N{u{x,y)) 



< M 



D 



N, 



H 



N'{u{x,y)) 



D 



<M^, 



oFi(l,\N{u{x,y)) 



< M 



D 



R, 



oFi(2,\N(Six,y))\)\N(Six,y))\ 



D 



<M'r. 



Using inequalities (15. 5p . (15. 6p together with notations (14. 7p . (14. Sp . (14.90 from equation (15.30 we 
can derive the estimate (Vx G [xi-i,Xi]) 



o u{x,yj^i] 



dx 



- 5Z / ^m|"(3^>?/)|^2/+ 



Vs-l 



(5.7) 



where 



i-i 



j/j-i 



+ X^ / N'-Jui,s\\uix,y)\dy+ / \Fix,y)\dy< 



Vs-l 



<h2j2Ni, 



|(fc)| 
u 



i-1 



T' l|(^)|l l|(°)| 



W, 



I l,S 



^2 5^Ar(.J|V||,J|M||_ + F||F|| < 



\l,S\\ \\l,S 



s=l 



s=l 



i-1 



(fc), 



(fc), 



< {Mn + M'^M^) /i2 J]|| « 11,^, + 1^11 F 



s=l 



u 



|W/ Ml II (^^11 II (^^ M 

h l^'?/) p. ' \\F\\ = \\F {x,y)\ 



(5.1 



IP,,,' 11^ II 11-^ v-'yy||D- 
Combining representation (15.40 with estimate (15.70 we obtain the inequality 

(k) (k) (k) 

ImII. . < (l + /ii/i2M') ||V||. , . + /ii/i2 M'M^M„ 11 V 11. . + 

I lUj — ^ '- ^ H-' W lU— Ij ' 1 ^ JV -n. "II lUj — 1 



i-1 



(fc), 



(fc), 



(fc), 



+/iiM^ <^ (M;v + M;M„) /i2 ^11 u II .^ + y|| f|| ^ + /ii/i2M^||f||. 



s=l 



l(fc)ll ll(fc)l 



(fc) 



ll(fcjll ll'-"'-'ll-l ^'^'' 
Denoting expression U^ • • LF by Ui.j, we can rewrite inequality (15. 8p in the form of 



(5.9) 



(fc) (fc) (fc) 

t/ij< (1 + /ii/i2M;j) Ui-i,j +hih2M'j^MRMu Uij-i + 
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+hiMR I /i2 (Mat + M'j^M^) 5J Uls +Y \ + h^h^MR. 

Using the method of mathematical induction it is easy to prove (see the proof of Theorem 
ED that 

(5.10) U,^j < fi,j, Vz G TJfi, J G MV2, 
where the real numbers fiij are defined by formulas (14.330 with 

(5.11) b = hMh) = hMR {h^M'^Mu + Y{Mn + M'^M^)) , 

a = l + hih2M'j^, c = hiMn{h2 + Y). 
Inequalities (IS.lOp together with Lemma [4.21 yields us the estimates 

(5.12) Ui,j < //ij < MRX{h2 + Y) exp((X + F)6i(/i2) + h^XM', 



Returning to the estimation of u{x, y) we get 
(5.13) 



E{h2) < E{1) = ai, Vz G 1, iVi, Vj e 1, A^2. 



{k)udefu(k\ ,|| ||(fc)|| ^ ||(fe)| 

n = M(x,?/j|g= max I "U 1 1 ^. ^. < 0"! 1 1 F I 

tg l,JVi 



i«,j 



Using estimate (I5.13P and equation (15.31) it is not hard to obtain the inequalities 



(5.14) 
where 



a u [x,y) 



dx 



D 



<Ya2\\F\ 



o u [x,y) 



dy 



D 



<Xa2\\F\\, 



(72 = ai{MN + M'j^Mu) + l. 
Combining inequalities (15.131) and (15.141) we get the following astimate 

(fc) . . M , M W, 



(5.15) 



u 



,,(fc) 11 iii/cj,, r , -j 

^ = II M (a;,y)||^^ < cr||F||, o" = max < ai, a2VX^ + Y'^> 



(k) 



Recalling the explicit formula for F{x,y) (see (15. 3p ) we can proceed with the estimation of 



I "U 11^ as follows 



(5.16) 



i(fc)i 
u 



k-1 



|(0)| 

u 



\(k-s),,.,,{s) 



U 



\ II "- ' II 

)lhll+ 



where 



fe-i 



(k-s-l) 



Er . /,. (0)/ N (k-s-l). .s 

^Ak_s-i[N]u±{x,y),..., u± {x,y))- 



s=0 



(0) 



(k-s-l) 



W. 



-Ak^s-i{N;u{x,y),..., u {x,y)) u{x,y) 



+ 



D 



+Ak [N; 11^11, ...,||u||j||n|| — ||'"||||'"||^(||'"||j), 



N{u) = J2 



uAu 



s=0 



To estimate the second term in the right side of inequality (15.160 we need the lemma stated 
below. 
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Lemma 5.1. 



As{N;uAix,y), . . . ,ujix,y)) - As{N; u{x,y),..., u{x,y)) 



< 



D 



<hAs{N,; 



|(0)| 

u 



u 



J, N,iu) = N'iu)u. 



Lemma [5. II is a partial case of Lemma 1 from [20] . 

Using Lemma 15.11 we can estimate the right side of inequahty fl5.16p in the following way: 

fe-i 



i(fc)i 
u 



,<a(j2'^''-^i^-^fA\i 



iC^-^)!! MlW 



U 



i)lhlli+ 



fc-i 



(5.17) 



s=0 



,(0) 



(fc-s-1),, ,,,(s). 



|(0)|| ||W|| ^~,,/||(0)| 



+ h2_^Ak^s~i[Ni] II u 11^,..., II u IIJII^IIi - II^IIiIpIIi^ (Ir^lli) )• 



s=0 



Let us consider the sequence of real numbers {t^fcjfclo defined by the formulas 

fc-i 



II (0) I 

vo=\\u\ 



Vk 



o- [y^ ^k-s {N;vo,..., Vk-s) Vs+ 



s=0 



fc-1 



(5.1^ 



+ ^Ak-s-i{Ni;vo, . . .,Vk-s-i)vs - VkN'{vo)voY k = 1,2,... 



s=0 



It is easy to see that for the sequence {vk}'kLo defined above the inequalities 

(k) 

(5.19) \\u\\^<Vkh'', A; = 0,1,... 

hold true. Assuming that the series 



(5.20) 



9{^) = E 



VkZ 



fc=0 



has a nonzero convergence radius, say R > 0, and g{R) < c>o, we immediately arrive at the 
inequality 

(5.21) VkR' < ^ 

with some positive parameters c and e. From inequality (15.211) it follows that the condition 

°2,||(fc)|| 
h < R is sufficient for the series 2^ Ki L to converge, that is, for the FD- method to converge. 

k=0 

Thus, to prove that inequality (15. ip holds for a parameter h, chosen sufficiently small, we have 
to investigate convergence of power series (I5.20p . 

Taking into account equalities (15.180 we arrive at the conclusion that function g{z) satisfies 
the nonlinear functional equation 

(5.22) (g^z)-vo){l + N'{vo)vo) = a[(^N{g{z)) - N{vo))g{z) + zN' {g{z)){g{z)) 

for all z G {—R, R). To prove that the radius of convergence of power series (I5.20p is nonzero, 
i.e, i? > 0, we have to consider the inverse function z = g~^. From equation (15.22P we can easily 
derive the explicit formula for z = z{g) : 



(5.23) 



49) 



{g - vo) (l + N'{vo)vo) - (^N{g) - N{vo))ga 



a 



N'{g)g^ 
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Taking into account that -2(^0) = we can easily find the value of z'{vq) : 

(5.24) z'ivo) = hm fM^fM = ^} 

g^vo g-vo aN'[vo)v^ 

Since function z{g) fl5.23p is holomorphic in some open neighborhood of the point g = Vq and 
z'{vo) > 0, we can conclude that there exists an inverse function z~^ = g which is holomorphic 
in some open interval {—R,R) (see [15]). Supposing that g{R) = 00 we get the contradiction 
(see (1E22])): 
(5.25) 

1 + N {vo)Vo = hm h a 

' g{z) 



N{g{z))-N{v^))+zN'{g{z))g{z 



+ CXD. 



This contradiction proves inequality (I5.2ip for some positive constants c and e which depend 
on N{u) and Vq only. Thus, we have that condition h < R provides the validity of inequalities 
dEH) and 

(5.26) 



i(fc)ii ^ c 



11 - A;i+^ \R^ 
Assume that h < R. Then inequality (15.11) allows us to consider the function 



00 , 
u(x. 



(k), 



y) = Y,u{x,y)eC{D). 



fc=0 



(fc). 



Furthermore, from equation (13.61) it follows that u{x,y) G C^'^(Pjj) and 



(5.27) 



d'^ u (x, y) 



dxdy 



< [Mn + M'j,Mu) 



u 



ik). 

f\ 



2Gl,iVi,jGl,A^2. 



Inequality (I5.27P together with (I5.26P imply that u{x,y) G C^'^(Pjj) and 



d'^ u{x,y) 
dxdy 



E 

k=0 



«2 ('^V N 

u[x,y) 
dxdy 



V(x,?/)gP,,„ iel,N,,j el,N2. 



The latter fact allows us to sum up equations (13.31) and (13.61) over k from 1 to 00 and this, 
taking into account the obvious equality 

00 k 

N{u{x, y)) = N{u{x, y)) u{x, v) = Y.Y. ^^-4^' 



(0) (k-s) (s) 

u{x,y),..., u {x,y)) u{x,y), 



k=0 s=0 



results in the equality 

(5.28) 



d'^ u{x,y) 00 



dxdy 



+ N{u{x,y)) = f{x,y), 



W{x,y) eDn{ix,y) \x^Xi,y^yj, ieO,NiJ g0,A^2}. 
Hence, we see that equality (I5.28P formally coincides with equation (12. ip . To obtain the identity 

u{x,y)=u{x,y), {x,y) e D 

00 00 

it is enough to remind that u(^0,y) = u{0,y) = (f){y), u(x,0) = u{x,0) = i/j{y) and to mention 
the fact that the solution u{x, y) to problem (12. ip . (12. 3 p is unique on D. 
Thereby, we have proved the following theorem. 
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Theorem 5.1. Suppose that the G ours at problem (12.11) . ( 12. 3 p satisfies the following conditions: 

oo 

(1) N{u) = uYj ^ku'', Vk G R, Vu e R; 

(2) ij{x) e CW (D^) n C {D^) , 0(y) G C(i) (D^) n C (D^) , /(x, y) G C(D). 

T/ien the FD-method described by formulas (13. ip . (13.20 . (13.30 . (13. 4p . (13. 6p . (13.70 converges su- 
perexponentially to the exact solution u{x,y) of the problem, i.e., the inequalities 

(5.29) \\u{x,yy ^(^^^)|li,z)< (^^i)L(^_/,) [r) ' ^ ^ N U {0} 

holds true, provided that h < R, where positive real constants c, R, e depend on the input data 
of problem (12. ip . (12. 3p only. 



6. Numerical example 
Let us consider the Goursat problem 

^^^ = e-(-), {x,y)eD, 
oxoy 

(6.1) u{x,0) = ^-ln{l + en, «(0, y) = | - ln(l + e^), 

where D = {^{x,y) | < x < 4, < y < 4}. Obviously, this problem satisfies the conditions of 
Theorem 15.11 It is not hard to verify that the exact solution to problem (16. ip is 

iz*(x,y) = ^^-ln(e" + e^). 



Using the FD-method described above we approximate the exact solution to problem (16.11) 

/ a finite sum fIXl 
Goursat problems: 



(fc) 
by a finite sum (13. ip with the terms u{x, y) satisfying the following recurrence system of linear 



92 u{x,y) l-exp(2 u{xi-i,yj-i)) (0) 

^^a^ + n ^ "^"' '^ ^ ' 

u[xi_i,yj_i) 

M(xi_i + 0,2/) =u{Xi-i -0,1/), m(x,7/j_i + 0) =m(x,|/j_i-0), 
(°)(^,0) = ^-ln(l + e-), ^n^(0,2y) = |-ln(l + e^), 

d^ u[x,y) l-exp(2 u (xi_i,|/j_i)) (fc) 
^^d^^ n ^ "^'''^^ 

-exp(2 n(a:i„i,2/j_i)) 2exp(2 ti(xi_i, 2/j_i)) \ (fc) (0) 

— ^^ ^ — + — ^^ — 1 ^(^-1' 2/i-i) <x^y)+ 

{u{xi-i,yj-i)) u{xi^i,yj^i) J 

(k) 
+ F{x,y), X G {xi^i,Xi), y G {yj-i,yj), 

f]f\ (k\ (k\ tk\ 

u{xi-i + 0,y) =u{xi-i -0,?/), u{x,yj-i + 0) =u{x,yj-i -0), ie l,A^i, j G 1,A^2 

(fc) (k) 

V(x, 0) = 0, V(0, ?/) = 0, A; = 1, 2, ... , 
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where function F{x, y) is defined by formula (15.31) and Xj, yj are defined according to (13. 2 p with 
X = Y = 4, Ni = N2 = 4; 20; 40; 80 (in terms of hi, /la we have that hi = h2 = 0.5; 0.2; 0.1; 0.05 
respectively). To monitor the accuracy of the method we use the function 

6{hi,h2,m) = II u{x,y,hi,h2) - u*{x,y)\\^. 
Table 1. The error of the FD-method as a function of the rank (m) and the mesh size {h) 





5(0.5, 0.5, m) 


5(0.2, 0.2, m) 


(5(0.1, 0.1, m) 


5(0.05, 0.05, m) 


m = 


1.0584498110834e-l 


1.3629587830264e-2 


4.3352923963359e-3 


1.0590305089182e-3 


m = 1 


2.0875237867244e-2 


5.6023714534399e-3 


2.0412391759766e-3 


3.6571985266298e-4 


m = 2 


1.5876742122176e-2 


1.7852756996399e-3 


3.1676334086428e-4 


4.3855534642367e-5 


m = 3 


8.7851563393853e-3 


1.7609349110392e-4 


2.0298330526525e-5 


1.5101132648798e-6 


m = A 


2.0883887112112e-3 


1.3084812991 115e-5 


1.1360343226132e-6 


4.6417353405381e-8 


m = 5 


2.9359063800745e-4 


5.1600423756071e-7 


5.6241341916952e-8 


2.7419476073821e-9 


m = Q 


1.7966735715635e-5 


4.4543002859311e-9 


4.2864933415766e-10 


5.3844093415473e-ll 


m = 7 


1.1298645650193e-6 


3.4087147338034e-10 


6.4824133982673e-ll 


3.7704350835172e-12 



The values of function 5(/ii, /i2, rn) presented in Table 1 show that the convergence rate of the 
FD-method increases as the mesh size decreases. Using the numerical data presented in the 
table it is not hard to verify that the function 6{hi, /12, m), as a function of m, (i.e., for the fixed 
values of parameters hi, /la) decreases exponentially as the rank m of the FD-method increases. 
This is in good agreement with Theorem 15.11 

7. Conclusions 

In the paper we have developed a numerical-analytic method for solving the Goursat problem 
for nonlinear Klein-Gordon equation. Under relatively general assumptions, we have proved 
that the method converges superexponentially, provided that the mesh size (h) is sufficiently 
small. From Theorem 15. II it follows that the accuracy of the FD-method can be increased either 
by increasing the rank m of the method or by decreasing the mesh size h. The latter conclusion 
has been confirmed by the results of the numerical example included in the paper. 

Though not discussed in the paper, the question of developing an efficient software implemen- 
tation for the proposed method is of great interest from both theoretical and practical point of 
view. We leave this question to a subsequent paper. 
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